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In  an  earlier  report  on  tidal  mixing  in  estuaries  [l]  an 
equation  for  the  diffusion  of  a solute  in  an  essentially  vertically 
homogeneous  estuary  in  which  the  major  influences  affecting  the 
movement  and  mixing  of  fluid  are  the  tides  and  river  flow  is 
derived,  namely 

o |S  + R 45  = D ±(x2  |S)  + q(x,t),  (1) 

at  9x  3x  Ox 


where  c is  the  time  average  of  the  estuary  cross-section  over 
a tidal  cycle,  assumed  here  to  be  independent  of  both 
x and  t; 

R is  the  river  discharge  and  has  dimensions  LJT"  ; 

2 -1 

D is  a positive  constant  having  dimensions  L T and 
is  defined  in  Eq.  116  of  reference  [1]; 


is  the  amount  of  solute  added  externally  and  has 
dimensions 


is  the  concentration  of  the  solute,  having  dimensions 
ML"3 ; 


x is  the  distance  along  the  estuary,  measured  downstream 
from  the  position  at  which  the  motion  in  the  estuary 
due  to  the  tides  is  assumed  to  be  zero* 

In  this  report  we  consider  several  solutions  of  (1);  in  the 

first  s is  the  concentration  of  salinity  in  the  estuary,  R is 

a function  of  time  and  q(x,t)  = 0;  in  those  following,  s is  the 


concentration  of  an  externally  introduced  solute,  R is  constant 
and  q is  assumed  to  depend  only  on  x. 

We  proceed  then  with  the  solution  of 


+ R - 

at 


(2) 


in  which  R = R(t),  and  D and  (a  + b)  are  constants.  It  is  assumed 
that  the  time  variation  of  R is  sufficiently  small  so  that 
(a  + b)  is  small  relative  to  the  other  terms  in  (2)  and 
hence  we  use  the  procedure  of  successive  approximations  in  the 


* Numbers  in  square  brackets  refer  to  the  bibliography  at  the 
end  of  the  paper. 
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The  boundary  conditions  to  be  placed  on  each  of  the  succes- 
sive approximations  of  s are  that  the  salinity  at  the  mouth  of 
the  estuary  he  d , a constant,  and  that  the  salinity  at  some 
point  upstream  be  zero.  Since  the  river  discharge  is  a function 
of  time  v»e  wish  to  permit  the  distance,  L,  from  this  point  to  the 
mouth  (the  effective  length  of  the  estuary)  to  vary  with  time, 
in  which  caso  it  is  convenient  to  choose  our  coordinate  system 
with  origin  at  the  mouth  of  the  estuary.  However,  with  this 
shift  in  the  origin  the  term  x^  in  (2),  which  is  proportional  to 
the  square  of  the  distance  moved  by  a particle  due  tc  the  tides, 

p 

must  be  replaced  by  (L  + x)  so  that  we  now  have 


= d£[cl  + x)2§f]  > 


(3) 


in  which  R and  L are  assumed  to  be  known  functions  of  time.  It 
must  also  be  assumed  that  R and  L do  not  vary  significantly  over 
a tidal  period  in  order  that  the  derivation  which  leafs  to  equa- 
tion (1)  remain  valid.  Physically,  this  appears  to  be  a reason- 
able assumption. 

Y/e  now  define  successive  approximations  sQ,  s^,  ...  to  tho 
solution  of  (3)  by 


D + x) 

O X 


(4  ) 


d ± kl  + x)2  ifn  ] » r + o 
9x  3x  9x 


dsv 


9s 


n=l 


at 


n= 1,2, ... 


(4.1) 


where  the  sn(x,t)  (n=0,l,2, , . .)  are  assumed  to  satisfy  tho 
boundary  conditions 

sn(0,t)  = b for  all  t,  (5) 


sn(-L(t),t)  = 


yJ 


for  all  t. 


(6) 
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Integrating  (4)  we  have 


D(L  + x)2  - Rs0  = f^t)  . (7) 

Here,  as  well  as  in  the  succeeding  equations,  the  f^Ct)  (1=1,2,.,.) 
will  be  used  to  denote  arbitrary  functions  of  time  which  will  be 
determined  by  the  boundary  conditions  (5)  end  (6),  Multiplying 
both  sides  of  (7)  by  the  integrating  factor 


B(L  + x)2 


_ R 

sdCl+xT 


and  integrating,  we  have 

s0(x,t)  = - £ih!  + f2(t)0^^  , 
R 

Imposing  (5)  and  (6)  on  this  solution  wo  have 


txW  = 0 , 


f2(t)  = o e1 


so  that 


R(1  _ _JL\ 
s0<x,t)-  = <>-eD  L Ux 


Integrating  (4.1)  we  have,  for  s1(x,t), 


D(L  + x)2  — 1 - Rs-,  = c \ --2.  dx  + f,(t)  . 

Ox  ^ at  d 


(10) 


The  lover  limit  in  the  integral  in  (10)  is  set  equal  to  -L  merely 
in  order  to  simplify  the  application  of  the  boundary  conditions. 
As  with  (7),  we  multiply  both  sides  of  (10)  by  the  integrating 


factor 


± — _ eJ 

D(L+x)2 


and  integrate,  obtaining 


s^(x,t)  a ce 


-R 

D(]j  *rX  ) 


— 1 — » e1 

D ( Xj+x  ) 2 


^ at 

-L 


dx  CX  - 


o 
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=L 


- -3—  + r At)*™**7  , 

H * 


(11) 


the  leer  limit  of  the  integration  again  being  chosen  in  order 
to  simplify  the  application  of  the  boundary  conditions.  Consider- 
ing the  double  integral  in  (11)  we  have,  upon  integrating  by  part^ 
x x 


D(L  + x)* 


“ e 


D(L+x) 


£!s 

at 


dx  dx 


iL 


1 

R 


R 

,D(L«) 


21 

at 


dx 


+ l 

+ R 


eDTL«T  |o  ax> 


-L 


(12) 


Substitution  of  (12)  in  (11)  then  gives 

x 


P -DdrKX  ) 

’l'*'*'"  ~ R - 


JL 


dx  - S. 


3s 

R j ST** 

-L 


8<1  - r±) 


+ 1 fiD  L L+x 
+ R 8 


3s, 


-R 


dx 


- + f .(t)oD^x\  (13) 


R 


-L 


From  (9)  virO  havo 


as-  _ fif/H,'  R'  , RL< 
“ D " t ■ - ' 


at 


L + * (L  + x)? 


(14) 


in  which  the  primes  denote  differentiation  with  respect  to  time, 

ds 

Thus  the  integrand  —-2  in  two  of  the  integrals  in  (13)  remains 

^ < < 

finite  over  tho  range  of  integration,  -L  = x *=  0,  Further, 
although  the  integrand  in  the  first  integral  on  the  right  hand 
side  of  (11)  becomes  infinite  at  x = -L, 

p x p 

Lim  oDTrSb  C^J  ^ dJC  „ n , 
r.-*  ~L  J 3t 


o 
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Thus  we  may  apply  boundary  condition  (6)  to  (13)  end  obtain 


f3(t)  = 0. 


Then  applying  (5)  to  (13)  we  have 


P 


fA(t)  = oeDL  , 


so  u.st 

S1(x , t)  = s0(x,t)  [1  + eDL 


rK  7 . B . 

eDnSn  tea  dJC 
St 


I(t4=  - 1) 


+ S-  I !fa  dx  - J.  a0'1*3'  ' L' 
or  at  or 


3 s 

dt2  to:1  • 


-L  -L 

We  consider  now  the  various  terms  in  (15).  From  (14) 


(15) 


=& 

JL  eDL 
OR 


e 


dTl+xT  8s 


X 

( 


= & [ (£>  - + ■ RLI  ■ 5 L •) 

DR  j L L + x (L  + x)2 


= £*(£)'  + .ffil  i0g 

nn'rV  + nw  -°b 


DR  L DR 


L + x 


(16) 


+ sL L(i  . — 3m-) 

D ' L L + x'  • 


Also  from  (14)  we  have 
R ° 


JSL  ©Di^  dx  - e 

or  j at  or  6 

-L 

R_  0 
= X ebL 


0 JCL«  - b (3=0 


31 


-L 


DR 


~R 


[*)  - Ji-  ♦ --R-LJ  leDU^J 

L L + x (l  + x)2 


-L 


P / R 

-A  ( [(a)«.  J! - + _au_).H?ST 

i-'H  } L u + x (L  + x)2 


dx 
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_ e_  [(B)'  . ax 

DR  L L + X 


ft  «.R  n ft  r _ft 

+ eiL  oBlLdTl^T  1 .4 1 

G j Q I • 


Noting  that  the  last  tv;c  terms  in  the  expression  directly  above 


cancel  each  other,  vje  have 


A ( dx 


B(_l_ 

c D L+x 
~ ah  e 


*’(  S* 


R 0 -R 

* reDL  ^ (1 a.)e^y  dx 

DIi ' 1/  Le  L+x 


R ? — R 

_ eD(L+xT  ( (1 a_)eDTE+xT  dx] 

L + x 


where 


u <«>' 

In  each  of  the  integrals  in  (17)  let  y = j so  that 

?>  R OD 

1 (1  - j^>ebT^57  dx  = | | (+  - ^o'ydy  . 

-L  tA,  V 

D(L+x) 


1 
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Integrating  the  first  term  in  the  integral  on  the  right  hand  side 
of  the  equation  directly  above  by  parts  we  have 


OD 


oo 


,-y  . 


ay  = 


e 


-y 


R 

dTl+xT 


co 

r 


e 


-y 


i — c'r 

j y ' * 


R R 

D(L+x)  D(L+x) 


so  that 
x 


-ft  -R 

(1  - _a_)Q^TL4xT  dx  = (L  + x)eD  ■L+x'  - (a  + |) 
L+x  D 


co 


e~y 


-L 


Substituting  (19)  in  (17)  we  have 


R 

DTb+xT 


o 

( 9s 


<jR 


i 


•rr2-  dx  - 
9t  rtR 


a<_i_  . i)  f 

D'L+x  l'  | 9s 


-L 


j ard*  = M>’j  - 


00 


+ (a  + §)[  eD^+x' 


e 


-y 


dy  - e 


R_ 

DL 


co 

-y 


e 


dy 


R 

nn£rr 


■ 


_R_ 

DL 


Finally,  substituting  (16)  and  (20)  in  (15) , we  have 


(20) 


3i(x,t)  = 3c(x't)ju  + wf  log  ife  | + ■ tfe5 

1 


+ li<1  + ^ 

* 


where 


oo 


,-y 


M(x)  = ex  - — dy  . 

y 


(21) 


(22) 


It  may  be  noted  that  values  of  the  function 

co 

- Ei(-x)  g I ^ dy 
J Y 

X 
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are  tabulated  (see  [2].) 
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For  the  case  in  which  the  mean  cross-sectional  area  c is 
a function  of  x the  formal  procedure  leading  to  s0  and  s^  is 
identical  with  that  just  given  although  the  details  of  the  inte- 
grations are  more  cumbersome.  The  term  D JL  (x2  «S.)  in  (1)  must 
° ox  9x 

then  bo  modified  slightly  and  may  be  obtained  from  Eq.  (97)  of 
reference  [ll  . 

Wo  now  consider  the  solution  of  (1)  for  the  case  in  which 
R is  constant  and  q(x,t)  » q(x)  so  that  we  may  assume  s(x,t)  a s(xX 
Then  (1)  becomes 

Dx2sh  + (2Dx  - R)s*  = - q(x)  , (23; 


in  which  primes  denote  differentiation  with  respect  to  x.  The 
chango  of  both  independent  and  dependent  variables  given  by 


z — — B— 
2Dx 


s(x)  S e“z4)(z>  , 
simplifies  (25)  to  the  particularly  nice  form 


,-20z 


V - <i>  q(^)  = • 


(24) 

(25) 


Us  in;;  the  method  of  variation  of  parameters  the  general  solution 

to  (25)  is  J{__ 

z 2DL 

'll  = '1>1  | 'l>2f(z)dz  + | \J>  xf(z)dz  4-  c-jC2  + o2°’"Z  > 

od  z 

where  and  c^  are  arbitrary  constants  and  andtj)^  ore  two 
solutions  to  the  homogeneous  equation  =0  such  that 

^1^2  “ ^ 2^1  = 1 » 
for  which  convenient  choices  are 
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Q 


i - e 


_ e 


~z 


1 


The  particular  choiae  for  the  limits  of  integration  in  (26)  are 
modo  in  order  to  simplify  the  application  of  the  boundary  condi- 
tions, Thus 


*<*>  - - a i 


R 

2DL 


:'2q(5D2)dz  - §T  j z"2e22q(~)dE  + c1°Z(26>1) 


jO 


+ c0e”Z  • 


Substituting  (24)  in  (26.1)  we  now  hove 
x L 


rS 


s(x)  =;  ^ 1 q(x)dx  + 1 eDx  eDxq(x)dx  + ^ + <?2 


JL 

,Dx, 


-R 

JDx 


(27) 


V/e  consider  two  special  coses,  namely 
(A.)  q(x)  — Qf>(x  - x0)  in  which  Q is  a constant,  0 < xD  < L, 

This  choice  for  q(x)  should  apply  to  the  introduction  of 
some  solute  at  one  point  (x  = xQ)  in  the  estuary.  Q is 
then  the  rate  of  discharge  of  this  solute  into  the  estuary 
(of  dimensions  MT"*). 

(D)  q(x)  = q in  which  q is  constant  over  the  range  0 < x < L, 
This  choice  for  q(x)  corresponds  to  introduction  of  solute 
uniformly  over  the  entiro  effective  length  of  the  estuary 
and  might  apply  in  the  investigation  of  ground  seepage.  The 
total  rate  of  discharge  of  the  solute  into  the  csturry  ovor 
the  range  0 < x < L is  then  qL, 

For  q(x)  as  given  in  (A),  (27)  becomes 

L 


x 

f 


s(x)  e £ 


6(x  - x )dx  + § o 
O R 


i r I 

Dx  1 „Dx  6(x  _ xQ)dx  + cx  + C2°DX 


o 


x 


f 
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Thus  if  we  imposo  the  boundary  conditions  s(0)  = s(L)  = 0 
(s(0)  being  set  equal  to  zero  since  there  is  no  mechanism  to 
carry  the  solute  upstream  of  the  point  x = 0,  at  which  point  the 
motion  in  the  estuary  duo  to  the  tides  is  zero;  and  s(L)  being 
set  equal  to  zero  because  of  the  cssontially  infinite  reservoir 
provided  by  the  ocean  at  the  mouth  of  the  estuary)  then 


ci  1 0 t 


"2  = - § . 

n 


so  that 


s(x)  = £ | Mx 

j 

o 


I 

xQ)dx  + £ oDx  i cDx  6 (x 


] 


S(i  . 1) 

x„)dx  -Sc01  * 

(28) 


which  may  be  written  in  the  form 
s(x 

a<i  . i) 

s(x)  = &[1  - eD  L x ] 


0 = x < x. 


(29) 


x0  < x 


< T 


For  q(x)  as  given  in  (B),  (27)  becomes 

=b  r j. 

s(x)  =flXaeDx  1 


R R 


e dx  + + c2 


Dx 


e 


Again,  imposing  tho  boundary  conditions  s(0)  = s(L)  = 0 we  have 

=a  ^ M. 

Dtt 

c^  * 0 since  lim  e 

x-to 


Dx 

e * dx  = 0 


and 


c2 


3L 

= . eDL 
R 
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L 
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S(JC)  = sf  - # e 


-i) 


Dx 


e 


R^ 

Dx 


dx 


(30) 


If  we  let  y = jijL  in  the  integral  in  (30)  then  we  have 


JL 

eDx  dx  s | 


_R. 

Dx 


R 

El 


y 

£-  dy 

y 


JL 

Dx 


= LeDL  . xeDx  + ft  Si  dy 


(31) 


R 

El 

after  integrating  by  parts.  Substituting  (31)  in  (30)  we  hove 

-R 


s(x)  = A e 


Dx 


eZdy 

y 


(32) 


_R_ 

DL 


It  may  be  interesting  to  compare  the  solute  distributions 
that  result  for  q(x)  as  given  in  (A)  and  (B)  respectively, 
assuming  that  the  rates  of  discharge  of  solute  throughout  the 
entire  effective  length  of  the  estuary  are  equal  for  the  two 
coses,  i.e.,  qL  =•  Q.  In  Fig.  1 Ra/Q  as  a function  of  x/L  is 
plotted  for  R/DL=0.8,  the  value  given  for  the  Raritan  River  in  a 
paper  by  Arons  and  Stommel  [ 3l • The  constant  which  they  denote 
by  F end  call  the  flushing  number  is  equal  to  R/DL  in  the 
notation  used  here.  Values  of  the  integral  appearing  in  (32)  are 
to  be  found  in  reference  [2], 
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Figure  1 presents  two  curves  for  the  case  in  which 
q(x)  =:  Q6(x  - x0)  j in  one  x0/L  = 0.4  and  in  the  other  x0/L  = 0.8, 
It  may  be  noted  in  particular  that  the  solute  density  is  appreci- 
ably greater  over  most  of  the  estuary  for  the  case  in  which 
x0/L  = 0.4,  It  is  to  be  hoped  that  this  sensitivity  of  the 
solute  density  upstream  of  the  external  source  to  the  position 
of  the  source  may  be  used  to  oheck  experimentally  the  validity 
of  the  solute  diffusion  equation  (1)  which  we  have  boon  using. 

From  (29)  we  may  note  that  the  maximum  value  of  s(x)  occurs 
at  x = xQ  in  the  case  in  which  q(x)  =■  Q6(x  - x0).  In  Fig.  2, 
Rs(x0)/Q  as  a function  of  x0/L  is  plotted.  A possible  application 
of  the  data  given  in  Fig.  2 is  the  determination  of  the  maximum 
distance  from  the  mouth  of  an  estuary  at  which  a known  source  of 
pollution  nay  bo  allowed  to  discharge  into  an  estuary  if  specified 
limits  aro  placed  on  the  permissible  density  of  pollution  in  the 
estuary. 
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fi|*  1.  fioluta  density  as  a 
function  of  distano* 

•Ions  astuary. 

A-*l  q(x)  - Q*(x-x0),  Xq/L  ■ 0.4 
A~2  qCx)  « QAtx-Xo),  x0/L  » 0.8 
B q(x)  ■ Q/L  0*x<L 


1.0 

ocean 
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rig.  2.  Mulaai  soluto  density  u 
• fonotion  of  position  of 
•xttrnsl  some*  In  ••txuorr 

qU)  * Q (x-*0) 


